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ABSTRACT 

Using two-dimensional hydrodynamical simulations, we investigate the interaction of a collimated fast 
wind (CFW) interacting with a spherical asymptotic giant branch (AGB) wind as the mechanism for 
shaping proto-planetary nebulae and young planetary nebulae. In particular, we compare our simulations 
to the observations of an evolved PPN with multiple, highly collimated lobes, CRL 618. We characterize 
our model CFW by three parameters: opening angle, velocity and mass-loss rate, and explore the 
dependence of the properties of the shell on the first two. For given opening angle and velocity, the 
mass-loss rate is chosen to give a shell velocity of about 150 km s _1 at the tip, similar to that seen in 
CRL 618. In our simulations, the shell dynamics is found to depend on the velocity of the fast wind: we 
, obtain a momentum-driven shell for a 300 km s^ 1 fast wind and a ballistic bow-shock driven shell for a 

^vq ■ 1000 km s _1 fast wind. The shell driven by the collimated fast wind is highly collimated, even though 

the AGB wind is spherical. Time variations in the velocity of the fast wind produce a series of internal 
7-H ■ shock pairs interacting with the inner surface of the shell. Due to radial expansion, the density of the 

^ ' internal shocks decreases with distance. 

Various emission diagnostics have been derived from our simulations. For a 300 km s _1 fast wind, 
the optical emission arises from both the shocked AGB wind and shocked fast wind, showing one or 
two bright bow- like structures at the tip of the lobe. However, for a 1000 km s _1 fast wind, since the 
shocked fast wind is much hotter, it emits mainly in X-ray emission; the optical emission forms only 
one bow-like structure at the tip associated with the shocked AGB wind. The position-velocity (PV) 
diagrams derived from our simulations all show a broad range of velocities at the tip. The detailed PV 
structure and velocity range at the tip depend on the shell dynamics and the relative contributions of 
i-^ ' the shocked fast wind and shocked AGB wind. 

We make a detailed comparison of our simulations to the observations of the relatively isolated north- 
western (Wl) lobe of CRL 618. We find that a 300 km s _1 collimated fast wind with an opening angle of 
10° can readily produce a highly collimated lobe similar to the Wl lobe, including the bow-like emission 
structure at its tip. However, our models have difficulty producing the bright emission structures seen 
along the body of the lobe. The [SII]A6716A/A6730A ratios at the tip of the lobe in all of our simulations 
are similar to that observed at the tip of the Wl lobe. The optical line ratios indicate a temperature 
• »"H . stratification in the tip, both for the simulations and observations, however, the temperatures at the tip 

of the lobe in our simulations are higher than observed. The position-velocity (PV) diagrams derived 
from our simulations are all qualitatively consistent with the current observations. The collimated fast 
wind in CRL 618 is unlikely to be steady and is not radiatively driven. 

Subject headings: stars: AGB and post- AGB — stars: winds and outflows. 



1. INTRODUCTION 

An investigation of the physical mechanism(s) respon- 
sible for the shaping of proto-planetary nebulae (PPNs) 
and young planetary nebulae (PNs) is crucial for under- 
standing the mass-loss processes during the end phases of 
the evolution of low- and intermediate-mass stars. Dur- 
ing the asymptotic giant branch (AGB) phase, these cool 
stars undergo intense mass-loss (with rate of up to 10 -4 
M & yr _1 ) via a dense, slow (~ 10 - 20 km s~ x ) wind (see 
e.g., review by Balick & Frank 2002), forming a circum- 
stellar envelope (CSE) of dust and gas around the star. 
During the PN phase, the central stars (now hot white 
dwarfs) lose mass via fast winds at speeds of ^1000-2000 
km s" 1 and rates of 10" 8 to 10" 6 M Q yr _1 (Balick & Frank 
2002). The fast wind catches up with the AGB wind and 
drives a shocked shell through it, forming a PN. This is the 
so-called interacting stellar winds (ISW) model commonly 



used to explain the formation of PNs (Kwok, Purton, & 
Fitzgerald 1978; Kahn 1983; Balick, Preston, & Icke 1987; 
Frank & Mellema 1994; Dwarkadas, Chevalier, & Blondin 
1996; Frank 1999). 

In the ISW model, the fast winds are spherical, and the 
AGB CSEs are assumed to have a density higher near the 
equator than near the pole, so that the shells driven by 
the fast winds expand fastest along the nebular symmetry 
axis. By varying the ratio of equatorial to polar density in 
the CSE, the ISW model is able to produce a variety of ax- 
isymmetric (e.g. elliptical or bipolar) shaped PNs (Balick, 
Preston, & Icke 1987; Frank & Mellema 1994). 

However, recent high spatial resolution and high 
dynamic-range imaging of many PPNs and young PNs 
with the Hubble Space Telescope (HST) shows that almost 
all of these have highly asphcrical shapes, with a signifi- 
cant fraction having highly collimated bipolar or multipo- 
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lar lobes (Sahai 2001; Sahai & Trauger 1998). Examples of 
objects with highly collimated lobes were also known from 
previous ground based observations (e.g., M 2-9, Ml-91, 
Ml-16, Corradi & Schwarz 1995; Schwarz et. al. 1997). 
Point symmetry, rather than axisymmetry, better charac- 
terizes the geometry of the majority of these objects [strik- 
ing examples are He 2-115 (Sahai & Trauger 1998), IC 4634 
(Sahai et al. 2003), and He 3-1475 (Borkowski, Blondin, & 
Harrington 1997)]. In view of such observations, Sahai & 
Trauger (1998) concluded that the ISW model could not 
explain the shaping of PNs, and proposed that collimated 
fast winds operating at the beginning of the post-AGB 
phase were the primary mechanism for initiating the shap- 
ing of PNs. 

Several nagging problems for the ISW model existed 
even before the HST imaging results became available. 
E.g., the shocked fast winds in the ISW model (with- 
out magnetic fields) did not converge to form stable, low- 
ionization microstructures (Icke, Balick, & Frank 1992; 
Frank & Mellema 1994; Dwarkadas & Balick 1998) like 
ansae, thin jets, or FLIERs (Balick et al. 1993), since at 
the high speeds of the fast winds, most windblown bubbles 
are "energy conserving" (Dyson & Williams 1980). MHD 
models with a toroidal magnetic field embedded in a nor- 
mal radiatively driven stellar wind are found to be able to 
produce such microstructures (Rozyczka & Franco 1996; 
Garcia-Segura 1997). However, the initial conditions in 
these models may not be appropriate (Gardiner & Frank 
2001). In addition, these models require radiatively driven 
winds which can not account for the very large momen- 
tum excesses observed in PPNs (Bujarrabal et al. 2001). 
Moreover, in order to produce bipolar PNs with collimated 
lobes and narrow waists, very high ratios of equatorial to 
polar densities of the AGB winds were needed (Balick, 
Preston, & Icke 1987; Mellema 1995). However, less than 
20% of the AGB envelopes show significant deviations from 
spherical symmetry (Bujarrabal & Alcolea 1991; Neri et 
al. 1998). Faint, roughly round halos can be seen around 
bright PPNs and young PNs, indicating that their progen- 
itor AGB envelopes are produced by isotropic mass-loss 
(Chu, Jacoby, & Arcndt 1987; Sahai 2001). The dense 
equatorial waist observed in many bipolar PPNs (e.g., Sa- 
hai et al. 1998; Kwok, Su, & Hrivnak 1998; Su et al. 1998) 
may actually result from the interaction between a colli- 
mated fast wind and the spherical AGB wind (Soker & 
Rappaport 2000). 

Thus, the evidence for the presence of CFW's in dying 
stars seems overwhelming, yet very little work has been 
done so far in order to investigate in detail how CFWs can 
shape the circumstcllar environment in these objects, and 
to constrain the physical properties of the CFWs through 
quantitative comparisons with observations. This paper 
makes a first attempt at addressing the above lack. We 
explore the interaction of a CFW with a spherical AGB 
wind using two-dimensional hydrodynamical simulations. 
The CFW is assumed to emanate radially from the source 
with a density decreasing from the pole to the equator. It 
is thus structurally different from a cylindrical jet, which 
emanates along the outflow axis with constant diameter 
and density. Cylindrical jet models are much more com- 
mon in the literature. These have been explored exten- 
sively in the case of YSOs (Blondin, Fryxell & Konigl 1990; 



Stone & Norman 1994; Biro & Raga 1994; Raga et al. 1995; 
Suttner et al. 1997), and to a much lesser extent, for PNs 
(Cliffe et al. 1995; Steffen & Lopez 1998). A CFW model 
has been used to produce the bipolar lobes associated with 
the luminous blue variable, r\ Car (Frank, Ryu, & David- 
son 1998). 

We support the results of our numerical simulations by 
comparison with analytical outflow models, and derive ob- 
servational diagnostics for comparison with observations. 
In §2 and 3, we describe the parameters and the numer- 
ical method used in the simulations. In §4, we present 
our simulations and compare them with analytical out- 
flow models. In §5, we compare our simulations with the 
observations of a young PPN, CRL 618. We summarize 
and conclude our work in §6. 

2. NUMERICAL SIMULATION PARAMETERS 

In our simulations, the AGB wind is assumed to be 
spherically symmetric with the density given by (in spher- 
ical coordinate) 

P ^ (!) 

where M a and v a are the mass-loss rate and velocity of the 
AGB wind, respectively. For our simulations to be com- 
pared with the observations of CRL 618, p a is calculated 
with M a = 3x 10~ 5 M Q yr" 1 and v a = 20 km s" 1 , the 
values comparable to those for CRL 618 (Knapp & Morris 
1985). After p a is calculated, however, v a is set to zero 
because it is very small compared to the velocity of the 
fast wind. This assumption of v a is appropriate in this 
paper, because we focus on the shell dynamics away from 
the equatorial regions. The temperature of the AGB wind 
is low and set to 10 K. 

At present, there is no consensus on how the fast wind 
is launched and collimated. One possibility is that the fast 
wind is launched from a magnetized accretion disk and star 
system by the magneto-centrifugal forces, as proposed for 
AGN jets (e.g., Blandford & Payne 1982) and protostellar 
winds/jets (e.g., Shu et al. 1994; Konigl & Pudritz 2000). 
These winds/jets are highly collimated by the magnetic 
fields. Keeping in mind our lack of understanding of the 
physical mechanisms driving the fast wind, we assume a 
simple parametrized description of the CFW density and 
velocity structure: 

P/ = P/o exp l-(8/0 f ) 2 } 

v^^expl-p;) 2 ]? (2) 
where pf and Vf a are the density and velocity at the pole 
(where = 0°), and 6f is the opening angle. The momen- 
tum per unit volume of the fast wind has a full width at 
half maximum of \/21n26>/ = 1.180/. For a given mass- 
loss rate of the fast wind (into two opposite sides of the 
central source), M/, the density at the pole is 

M f 

Pfo - 



Amr^VfoK 



(3) 



where rj is the launching radius of the fast wind and 

r /2 o 

K= / cxp[-2(— ) 2 ] sin (9 d<9 (4) 
Jo Of 

In our simulations, p/ D is calculated from Mf with 77 as- 
sumed to be 5x 10 15 cm (333 AU), so that our fast wind is 
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characterized by three parameters: Mf, Vf a and Of. Orig- 
inating near or from the post-AGB star, the fast wind is 
assumed to have a temperature of 10 4 K. 

The parameters of the fast wind in our simulations are 
listed in Table 1. Observationally, the fast wind is found 
to have a velocity up to <~ 2000 km s _1 (e.g., He 3-1475, 
Sanchez Contreras & Sahai 2001). Here, Vf a is set to two 
velocities, 300 and 1000 km s _1 . In the simulation with 
Vf — 1000 km s _1 , Vf is allowed to increase linearly from 
to 1000 km s _1 in the first 10 years to avoid a very strong 
shock at the beginning. As we will see later, these two ve- 
locities of v f result in qualitatively different dynamics in 
the simulations. Of is set to 10° and 20° in order for our 
simulations to be compared with highly collimated PPNs 
and young PNs. Mf ranges from 2.5xl0~ 7 to 10~ 5 M Q 
yr^ 1 . It is chosen to give a shell velocity at the tip of about 
150 km s _1 , which is comparable to that seen in CRL 618 
(Sanchez Contreras, Sahai, & Gil De Paz 2002, hereafter 
SCSG02). 

The fast wind may have periodical variations in density 
and velocity, as suggested in e.g., CRL 618 (SCSG02). 
Such variations could result from variations in the accre- 
tion process as suggested for protostellar winds (Hartmann 
& Kenyon 1996). In this case, the density and velocity are 
assumed to vary in such a way that the mass-loss rate of 
the fast wind is constant over the variation. Hence, the 
fast wind is assumed to have 

P/=P//(1 + Asin— ) 

v , /=v / (l + ^sin^) (5) 

where A and P are the amplitude and period of the vari- 
ation. Such variations in density and velocity have also 
been used in simulations of protostellar winds /jets in or- 
der to produce knotty and bow-like structures along the 
outflow axis (see e.g., Biro & Raga 1994; Lee et al. 2001). 

3. EQUATIONS AND NUMERICAL METHOD 

The two-dimensional hydrodynamic code, ZEUS 2D, is 
used to solve the equations of hydrodynamics, 

! + V.(pv)=0 

^ + v • h = - Vp + / 

de 

— + V • (ev) = -pV • v - n 2 A (6) 

where p, v, p, e, and n are the mass density, velocity, 
thermal pressure, internal energy density, and hydrogen 
nuclei number density, respectively An ideal gas equation 
of state p = (7 — l)e with 7 = 5/3 is used for the ther- 
mal pressure. Helium is included as a neutral component 
with n(He) = 0.1 • n, so that n = p/(1.4 • ma), where 
rriH is the mass of atomic hydrogen. In our simulations, 
the undisturbed AGB wind is assumed to be static. As 
a result, a small fixed force per unit volume, /, is intro- 
duced in our simulations to hold the undisturbed AGB 
wind from thermal expansion. Here / = — ^2., where p a 
is the thermal pressure of the undisturbed AGB wind. A 
is the optically thin radiative cooling function for inter- 
stellar gas, with the function at high temperature from 



MacDonald & Bailey (1981) and the function at low tem- 
perature from Dalgarno & McCray (1972). We assume an 
equilibrium cooling function, in which the ionization frac- 
tion of hydrogen is calculated by equating the ionization 
rate with the recombination rate of hydrogen. As a result, 
the cooling function decreases rapidly at about 10 4 K as 
the ionization fraction of hydrogen decreases. Molecular 
cooling is not included. 

Scalar tracers, c, c a , and c w , are also included in the 
simulations to track the fast and AGB winds, c tracks the 
fast wind, it is one for the fast wind, zero for the AGB 
wind, and a value between one and zero for a mixture of 
the fast wind and AGB wind. c a and c w trace respec- 
tively the flows of the AGB wind and the fast wind in the 
strongly interacting regions near the pole. Thus, c a (c w ) 
is one for the AGB (fast) wind emanating within the an- 
gular segment 5°< <7.5°, a value between one and zero 
for a mixture of this material with the rest, and zero for 
the rest. 

Simulations arc performed in spherical coordinates with 
a computational domain of dimensions (r, 0) = (5 x 
10 15 to 10 17 cm, to b ) and a uniform grid of 400 x 400 
zones. This gives a resolution of 2.37 x 10 cm in r. Of, is 
a value between 0.5 and ir/2 depending on the angular ex- 
tent of the simulations, giving a resolution of 1.25 x 10~ 3 to 
4xl0~ 3 radians in 0. Reflecting boundary conditions are 
used along the inner (where — 0°) and outer (where 
= Of,) boundaries, while outflow boundary conditions are 
used along the outer r boundary. Inflow boundary condi- 
tions are used for inner r boundary to introduce the fast 
wind into the AGB wind. The results of the simulations 
are presented in cylindrical coordinates (z,R), with the z- 
axis being the outflow axis. 

4. RESULTS 

In this section, we present the results of our simulations 
as the fast wind propagates to a distance of about 10 17 cm, 
a length similar to that of the lobes in CRL 618. Observa- 
tionally, a significant amount of optical emission is believed 
to come from a shell surrounding the fast wind, therefore, 
we concentrate on the shell dynamics in the simulations. 
Model 1 is considered as a standard model because it can 
produce a shell with a collimation similar to that of CRL 
618. We discuss the effects of changing the opening angle 
and velocity of the fast wind on the structure and kine- 
matics of the shell. We use analytical models to support 
the results of our simulations. At the end of this section, 
we discuss how the velocity of the fast wind changes the 
shell dynamics in the simulations. 

4.1. Model 1: Standard model 

The simulation of model 1 is presented in Figure 1. The 
fast wind is 177 years old. In this model, Mf = 2.5 x 10~ 6 
Mq yr" 1 , Vf = 300 km s"\ and f = 10°. As the fast 
wind blows into the AGB wind, a shell is formed around 
the fast wind. The shell is highly collimated with a length 
to width ratio of about 6. It consists of shocked AGB 
wind (forward shock) and shocked fast wind (backward 
shock) separated by a contact discontinuity (mixture), as 
indicated by the c contours. Since the heat generated by 
the fast wind's kinetic energy is efficiently converted into 
escaping radiation, the shocked fast wind becomes a thin 
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inner shell accelerating the shocked AGB wind in the outer 
shell. As a result, the shell is a "momentum-conserving" 
wind-driven shell. The motions of the shocked fast wind 
and shocked AGB wind are different. The shocked fast 
wind moves along the shell toward the pole, producing a 
focusing effect that in turn produces a high-density jet- 
like structure at the tip, as indicated by the c w contours 
(see Canto, Tenorio-Tagle, & Rozyczka 1988; Frank, Bal- 
ick, & Livio 1996b). The detailed structure of the jet-like 
structure, however, depends on the numerical resolution. 
The shocked AGB wind moves along the shell toward the 
equator, as indicated by the c a contours. This motion of 
the shocked AGB wind depends on the relative velocity 
between the shell and undisturbed AGB wind and thus 
will still be important even if the undisturbed AGB wind 
has a radial expansion of 20 km s _1 . However, this motion 
is not able to produce mass accumulation at the equator 
due to high density (thus inertia) of the shell material near 
the equator. In the simulation, no Kelvin- Helmholtz insta- 
bility is seen along the contact discontinuity because the 
velocity changes smoothly within the shell and the resolu- 
tion in our simulation is low. The velocity decreases from 
the shocked fast wind in the inner shell to the shocked 
AGB wind in the outer shell, indicating that the momen- 
tum of the shocked fast wind is not transferred immedi- 
ately to the shocked AGB wind due to the thermal pres- 
sure in the shell. The shell velocity, v s , decreases with 
latitude from the pole to the equator but increases with 
time. The velocity at the tip where 9 = 0° is found to be 
v so = 138(1 + 0.00165i yr ) km s _1 , where t yr is the outflow 
age of the fast wind in units of years. At the age of 177 
years, v so is 178 km s _1 . 

The fast wind itself cools rapidly as it propagates away 
from the source due to adiabatic expansion and radiative 
cooling. For the temperature of the fast wind, Tf, we have 
dT f 4 Tf 2nA 
dr 



(7) 



3 r 3kvf 

Due to radial expansion, the number density n = 

lArriH ^ r ) Thus, 

dT L = _AT L 
dr 



3 r 



Pf A 2 
2.lniHkvf r 



(8) 



Integrating this equation numerically, we find that the fast 
wind cools down from 10 4 to 10 2 K after traveling a dis- 
tance of 1.2x 10 16 cm from the source, consistent with the 
simulation. The majority of the shell has a temperature 
of about 10 4 K at which the cooling function decreases 
rapidly (see §3). The temperature of the shell is highest 
at the shell boundaries near the pole where the material is 
recently and strongly shocked, about 10 5 K for the shocked 
AGB wind and a few 10 4 K for the shocked fast wind. The 
temperature decreases away from the pole as the shell ve- 
locity decreases. The temperature also depends on the 
radiative cooling and thus the density. Near the source 
where the density is high, the shocked AGB wind cools 
to less than 100 K (due to the cooling of atomic gas at 
low temperature in our cooling function, see Dalgarno & 
McCray 1972). The dense jet-like structure at the tip also 
cools to a very low temperature compared to its surround- 
ing shell material. 

In order to further investigate the shell dynamics, we 
present the transverse velocity, vr, and longitudinal ve- 
locity, v z , of the shell material in Figure 2. The velocity 



structures of the newly shocked fast wind and the newly 
shocked AGB wind arc different. The transverse velocity 
of the newly shocked fast wind decreases with increasing 
z. It becomes negative near the pole but drops rapidly to 
zero at the tip, indicating that the newly shocked fast wind 
is focused toward the pole and shocked strongly at the tip. 
At the tip, the sound speed of the newly shocked fast wind, 
c s , is 12 km s~ 4 . According to Bernoulli's theorem, the 
transverse velocity of the newly shocked fast wind is lim- 
ited to V3cs or 20 km s _1 for 7 = 5/3. On the other hand, 
as we go from the tip to the central source (decreasing z), 
the transverse velocity of the newly shocked AGB wind 
first shows a rapid increase in the tip (indicating that the 
newly shocked AGB wind at the tip is squirted sideways 
along the shell) and then shows a steady linear decrease. 
Since the momentum of the fast wind is mainly along the 
z-axis, the longitudinal velocity of the newly shocked AGB 
wind is very low except near the tip. 

In the following, we compare the shell structure and 
kinematics in the simulation of model 1 to a simple ana- 
lytical momentum-driven shell model given by Lee ct al. 
(2001), which was found to be able to describe the shell 
dynamics in similar but isothermal wind simulations. In 
this analytical model, the thermal pressure in the shell and 
the motion of the shocked fast and AGB winds along the 
shell are not considered. In this case, the ram pressure of 
the fast wind is exactly balanced by that of the AGB wind, 
so that the shell expands outward with a constant velocity 
along any radial direction. Hence, the shell velocity (see 
Figure 3) is, 



i + n 



and the shell radius is 



r s {6) =v s (6)t + r f 



(10) 



where 77 



1 M 



-TT- — - Therefore, the shell velocity 



= r f ps 

r 2 p a K M a Vf ' 

at the tip, v so , is expected to be 138 km s _1 , in agreement 
with the initial shell velocity as the fast wind first impacts 
the AGB wind in the simulation. However, in the simula- 
tion, the shell velocity increases with time as the shocked 
fast and AGB winds flow along the shell. As a result, the 
predicted shell structure (the dashed line in Figure la) is 
smaller than that in the simulation. 

The velocity of the newly shocked fast wind can be de- 
rived from the shell velocity. However, since the shell ve- 
locity increases with time in the simulation, a mean (av- 
eraged over time) shell velocity, v s , is used to derive the 
velocity of the newly shocked fast wind. The mean shell 
velocity is v s (9) = (r s (#) — r/)/i. As expected, its R and z 
components (the solid lines in Figure 2) lie between those 
of the newly shocked fast wind and the newly shocked 
AGB wind. The velocity of the newly shocked fast wind 
is 

v s/ =v/-^^ — cosa[sin(a+6>)R+cos(a+fl)z] (11) 
7 + 1 

where a is the angle between the shell normal and the 
mean shell velocity (see Figure 3) - note that the depen- 
dence of various parameters on 9 is not explicitly shown in 
the equation. Its R and z components (the dashed lines in 
Figure 2) are found to match reasonably well those of the 
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newly shocked fast wind in the simulation, except near the 
tip where the shocked fast wind shocks with itself, produc- 
ing a thermal pressure gradient that reduces the transverse 
velocity near the tip. With this equation, the transverse 
velocity of the newly shocked fast wind is found to become 
negative as 9 < 6 C , where 

sin6> c = 2( ~7 - ^ cosasinfa + 6 C ) (12) 
w/(7 + l) 

As a result, the shocked fast wind with 6 < 9 C or 3.4° is 
focused toward the pole, producing a focusing effect that 
produces a shell longer than the momentum-driven shell 
in the simple analytical model. 

4.2. Model 2: Effects of changing the opening angle 

The effects of the opening angle of the fast wind can be 
studied by comparing the simulation of model 2 to that of 
model 1. The simulation of model 2 is presented in Figure 
4. The age is 188 years. In this model, M f = 1 x 1CT 5 
M yr" 1 , v fo = 300 km s _1 , and 6 f = 20°. The mass-loss 
rate is adjusted to give the same pf and thus the same 
initial shell velocity at the tip as model 1 (sec Equation 
9). The shell dynamics in this model is qualitatively the 
same as that of model 1. However, driven by a fast wind 
with twice opening angle, the shell has a length to width 
ratio of about 3, half of that in model 1. Hence, under 
our assumption of the fast wind, the collimation (length to 
width ratio) of the shell increases roughly linearly with de- 
creasing opening angle of the fast wind. In addition, since 
the shell velocity decreases slower away from the pole, the 
temperature at the shell boundaries also decreases slower 
away from the pole. Moreover, since the angle between the 
shell normal and the shell velocity is everywhere smaller, 
the focusing effect of the shocked fast wind is less. As a re- 
sult, the shell velocity increases with time slower than that 
of model 1 , so that the shell structure can be described by 
the simple momentum-driven shell model (see the dashed 
line in Figure 4a). The shell velocity at the tip is found to 
be v so = 138(1 + 0.0004V) km s"\ or 148 km s" 1 at the 
age of 188 years. 

4.3. Model 3: Effects of changing the velocity 

The effects of the velocity of the fast wind can be studied 
by comparing the simulation of model 3 to that of model 
1. The simulation of model 3 is presented in Figure 5. 
The age is 180 years. In this model, M f = 2.5 x 10~ 7 M Q 
yr^ 1 , Vf = 1000 km s" 1 , and Of = 10°. The mass-loss 
rate is chosen to give an initial shell velocity of about 150 
km s _1 at the tip, similar to that of model 1 (see Equa- 
tion 9). Although the opening angle of the fast wind is 
the same, the shell is broader than that of model 1. Un- 
like model 1, the shell is mainly composed of the shocked 
AGB wind, as indicated by the c contours. The shocked 
fast wind does not cool fast enough and expands into a 
low-density but hot, high-pressure and turbulent cocoon 
surrounding the unshocked fast wind, as indicated by the 
c w contours. Being compressed by this high-pressure co- 
coon, the fast wind itself becomes a ripple jet surrounded 
by a dense sheath of compressed fast wind, similar to that 
seen by e.g., Frank & Mellema (1996a). The jet interacts 
with the AGB wind, producing a rippled shell with a con- 
vex curvature at the tip. The shocked AGB wind in the 



convex curvature moves toward the symmetry axis. How- 
ever, since the shocked AGB wind cannot flow across the 
symmetry axis in the simulation, some of it is deflected 
in a backward direction toward the source, producing a 
linear and dense structure along the axis inside the shell, 
as seen in other simulations (Blondin & Lundqvist 1993; 
Dwarkadas & Balick 1998). As a result, this linear and 
dense structure is likely an artifact produced by the nu- 
merical effect of the simulation. Since the fast wind itself 
is compressed by the shocked fast wind and bent toward 
the z-axis, the density of the fast wind does not decrease 
with distance as fast as the density of the AGB wind, so 
that the shell velocity at the tip increases with time. At 
the tip, the shell velocity is about 220 km s -1 and the 
temperature of the shell is about 10 6 K, much higher than 
that in model 1. On the other hand, the temperature of 
the shocked fast wind is a fewxlO 7 K at the head of the 
fast wind and decreases to a fewxlO 6 K as the shocked 
fast wind expands into the cocoon. Since the shocked fast 
wind expands sideways at the head, there is no focusing 
effect in this simulation. 

The transverse velocity, vr, and longitudinal velocity, 
v z , of the shell material are presented in Figure 6. Since 
we are interested in the shell dynamics, the shocked fast 
wind in the cocoon is excluded in this figure. The negative 
vr is associated with the convex curvature at the tip. To 
avoid confusion, the convex curvature is not included in 
the figure of v z . Unlike model 1, the transverse and longi- 
tudinal velocities both show a convex spur structure with 
the highest velocity at the tip, indicating that the shocked 
AGB wind at the tip, due to its high thermal pressure 
and the thermal expansion of the shocked fast wind at the 
head, is squirted sideways along the shell. This is similar 
to that seen in the simulations of a cylindrical jet (see e.g., 
Lee et al. 2001). 

In the following, we compare the shell structure and 
kinematics in the simulation of model 3 to a simple ana- 
lytical ballistic bow shock model in Ostrikcr et al. (2001), 
which was found to be able to describe the shell dynam- 
ics in the simulations of a cylindrical jet. In the ballis- 
tic bow shock model, the shocked AGB wind in the shell 
is assumed to be accelerated by an initial impulse at the 
head - the thermal pressure forces of the shocked fast wind 
and shocked AGB wind are ignored subsequent to the ini- 
tial impulse. At the tip, the shell velocity, v so , is about 
220 km s^ 1 and the sound speed is about 65 km s _1 in 
the simulation. Therefore, vr is limited to the smaller 
value of VSc s and v so /2, or 110 km s _1 , consistent with 
the ballistic bow shock model. Using Equations 18-22 in 
Ostriker et al. (2001) and assuming a mean jet radius of 
Rj = 4 x 10 15 cm for the ripple jet (see Figure 5) and 
= 4 (a typical value used for the initial impulse at the 
head, see Lee et al. 2001; Ostriker et al. 2001), we calcu- 
lated the shape of the shell (the dashed line in Figure 5), 
the mean shell velocity averaged over the shell thickness 
(the solid lines in Figure 6), and the immediate postshock 
velocity of the AGB wind (the dashed lines in Figure 6). 
As can be seen, the shape and the transverse velocity of 
the shell in the simulation can be reasonably reproduced 
with the ballistic bow shock model. However, the simple 
ballistic bow shock model cannot reproduce the longitudi- 
nal velocity that well. In the simulation, the material in 
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the shell consists of two components, one with a velocity 
close to the immediate postshock velocity and the other 
with a velocity close to the mean shell velocity. 

4.4. Model 4: Effects of the time variations 

The effects of the time variation of the fast wind can be 
studied by comparing the simulation of model 4 to that of 
model 1. The simulation of model 4 is presented in Figure 
7. The age is 166 years. In this model, the fast wind has 
the same parameters as that in model 1 but has a time 
variation in the density and velocity with A — 0.5 and 
P = 22 years (see Equation 5). The structure and veloc- 
ity of the shell in this model are similar to those of model 
1. There is also a jet-like structure at the tip. Inside the 
shell, however, there are periodic internal shock pairs (i.e. 
forward and backward shocks), which are clearly seen in 
the temperature distribution, propagating away from the 
source at the velocity vj. The wings of the internal shocks 
grow with time and eventually impact the inner surface 
of the shell. The curvature of the internal shocks are de- 
termined by the angular distribution of the fast wind's 
velocity and thus well described by (dashed lines in Fig. 
7a) 

r(6)=v f (0)t t + r f (13) 

where U is the age of each internal shock. The thickness 
of the internal shock pairs (i.e the separation between the 
forward and backward shocks) increases with the distance 
due to thermal expansion. The thickness along the z-axis 
(6 = 0°) is found to be 



giving a kinetic energy injection rate along the z-axis of 
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where c s = 11 km s , the sound speed at the shock 
boundaries in the simulation. The density of the internal 
shocks decreases with the distance due to radial expan- 
sion. The shocked fast wind also moves along the internal 
shock surfaces, as indicated by the c w contours. The lead- 
ing and the internal shocks all show highest temperature 
at the boundaries where the material is newly shocked. 

4.5. Shell dynamics 

It is clear that the velocity of the fast wind is a crucial 
parameter in determining the shell dynamics, e.g., whether 
the shell is momentum- driven or driven by thermal pres- 
sure. This was recognized by Dwarkadas & Balick (1998) 
who carried out simulations of an ISW model in which 
the fast wind velocity ramped up from 25 to 1000 km 
s" 1 - their simulated nebulae evolve through an initial 
momentum-conserving phase before entering the energy- 
conserving stage. In our simulations, the parameters are 
chosen to give a shell velocity of about 150 km s -1 at the 
tip, comparable to that seen in CRL 618. For a given 
shell velocity at the tip, v so , the velocities of the shocked 
AGB wind and shocked fast wind at the tip are v so and 
Vf — v so = v s f , respectively. In order to balance the ram 
pressures of the shocked fast wind and shocked AGB wind 
at the tip, the density of the fast wind along the z-axis 
(8 = 0°) is required to be 
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For a 300 km s _1 fast wind, the shocked AGB wind and 
shocked fast wind both have a velocity of about 150 km 
s _1 at the tip. The cooling is fast enough to remove the 
heat obtained from the fast wind, so that the shell is thin 
and momentum-driven. For a 1000 km s _1 fast wind, how- 
ever, the shocked fast wind has a velocity of 850 km s _1 
at the tip. Thus, the fast wind is required to be more 
tenuous than that of 300 km s _1 fast wind (see Eq. 15). 
However, the energy injection rate is about the same (see 
Eq. 16), so that the shocked fast wind does not cool fast 
enough and expands into a hot cocoon. The fast wind it- 
self becomes a ripple jet surrounded by a dense sheath of 
compressed fast wind. As a result, the shell is essentially 
driven by a jet. At a velocity of about 150 km s , the 
shell is radiative and its kinematics can be described by a 
ballistic bow shock model. 

4.6. Comparison with previous analytical models 

Recently, Soker (2002) has also developed analytical 
models of highly CFWs in order to reproduce highly colli- 
mated shells in PPNs and PNs. He defines two cases: (i) 
the "slow-propagating" jet, in which the shell at the tip 
proceeds at a speed v so <C Vf , which can be compared to 
our simulations of model 3; and (ii) the "fast-propagating" 
jet, in which v so ~ Vf a . In case (ii), the density of the gas 
in the jet is much larger than that of the AGB wind; this 
condition is not met in any of our simulations. In Soker's 
case (i) above, the shocked fast wind is assumed to be 
adiabatic and expands sideways behind the shell of the 
shocked AGB wind, similar to that seen in our simulation 
of model 3. He also assumes that the shell is continuously 
accelerated by the thermal pressure force of the shocked 
fast wind and expands laterally perpendicular to the out- 
flow axis. In our simulation, however, the ram pressure 
in the shell is much greater than the thermal pressure of 
the shocked fast wind except at the head. As a result, 
the force due to the thermal pressure of the shocked fast 
wind can be ignored subsequent to an initial impulse at the 
head, so that the shell dynamics is in the ballistic limit. 
In addition, unlike his model, the shell in our simulation 
also expands significantly in the z direction (see Figure 6) 
as the momentum of the fast wind is transferred to the 
shocked AGB wind. As a result, his model predicts a lobe 
much wider than that in our simulation (model 3). 

5. COMPARISON WITH OBSERVATIONS 

In the following, we compare our simulations to the ob- 
servations of CRL 618. CRL 618 is a young PPN located 
at a distance of 0.9 kpc (Goodrich 1991). It shows sev- 
eral narrow lobes at different orientations in HST images 
(Trammcll 2000), perhaps resulting from multiple ejec- 
tions at different orientations. The general structures of 
the different lobes are similar. We compare our simulations 
to the northwestern (Wl) lobe, which is separated from 
other lobes and shows clear position-velocity diagrams in 
the groundbased observations. This lobe has an inclina- 
tion of about 30° into the plane of the sky (SCSG02) and 
a length of about 10 17 cm. Since molecular cooling is not 
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included in our simulations, we only compare the optical 
emission. In the simulations, the emission is assumed to be 
optically thin and arising from gas in local thermal equi- 
librium (LTE) (see Appendix). We assume a distance of 1 
kpc and an inclination of 30°, the values similar to those 
of CRL 618. 

5.1. Morphology 

The Wl lobe of CRL 618 is highly collimated with a 
length to width ratio of about 7. In optical emission, sev- 
eral bow-like (or ripple- like) structures are seen within the 
body of the lobe (see Figure 8). The [SII]A6716A/A6730A 
ratio is about 0.5 and roughly constant in the lobe 
(SCSG02). We compare our simulations to these features 
in order to determine which model can better produce the 
morphology of the lobe. 

With a highly collimated fast wind, model 1 is able to 
produce a narrow lobe similar to the Wl lobe. In this 
model, the velocity of the shocked AGB wind is higher 
than that of the shocked fast wind, so that the opti- 
cal emission arises mostly from the newly and strongly 
shocked AGB wind, showing a single bow-like structure at 
the tip (see Fig. 9). From [OIII] to [Nil] to [SII] to [01], 
the emission traces the material at lower and lower tem- 
perature, showing more and more extended bow-like emis- 
sion structure. At the tip, the electron density is about 
10 4 cur 3 and the [SII]A6716A/A6730A ratio is about 0.5, 
similar to that seen in the Wl lobe. The emission from 
the newly shocked fast wind could be as strong as that 
from the newly shocked AGB wind, showing two bow-like 
structures at the tip, as that seen in model 2 (see Fig. 
10), in which the velocities of the shocked AGB wind and 
shocked fast wind happen to be comparable, about 150 
km s -1 . Due to its high density, the material between the 
two bow-like structures cools down quickly and radiates 
very weakly in optical emission. As a result, a 300 km s _1 
fast wind with a small opening angle can produce one or 
two close bow-like structures at the tip similar to that seen 
in the Wl lobe, however, it has difficulty producing other 
bow-like structures in the Wl lobe. 

On the other hand, a model with a 1000 km s _1 fast 
wind has difficulty producing a narrow lobe similar to that 
of the Wl lobe. Although it has the same opening angle 
of the fast wind as that in model 1, model 3 produces a 
lobe much broader than the Wl lobe. We have performed 
a simulation using a fast wind with an opening angle half 
of that in model 3 but with the same velocity and den- 
sity along the axis, however, the lobe in that simulation 
is only 15 percent smaller than that of model 3, and thus 
still much broader than the Wl lobe. As a result, due 
to thermal expansion of the shocked fast wind, the col- 
limation of the lobe produced by a high speed fast wind 
increases slowly with decreasing opening angle of the fast 
wind. In model 3, the optical emission shows a bright 
bow- like structure at the tip (see Fig. 11). The emission 
becomes weak away from the tip. Except near the tip, 
the [OIII] emission is from the inner shell while the [OI] 
emission is from the outer shell. The [SII]A6716A/A673oA 
ratio is everywhere about 0.5, similar to that seen in the 
Wl lobe. There are ring-like emission structures arising 
from the ripple in the shell in the simulation (see Fig- 
ure 5). However, these ring-like structures are much more 



rounded than the observed bow-like structures in the Wl 
lobe, suggesting that the bow-like structures are not asso- 
ciated with the ripple in the shell. There is also a cone- 
like structure along the axis associated with the dense 
sheath of compressed fast wind around the jet (see Fig. 
5). Since this structure is not seen in the Wl lobe, the hot 
and high-pressure cocoon that produces the dense sheath 
of compressed fast wind is probably not present in the 
Wl lobe. As a result, we conclude that the velocity of 
the shocked fast wind in the Wl lobe is probably smaller 
than that in model 3. Since the hot and high-pressure 
cocoon radiates in X-ray emission, future X-ray observa- 
tions will be able to check the existence of the hot cocoon 
in the Wl lobe. For future comparison, we present the 
X-ray emission and spectrum in the 0.2 -1.5 keV energy 
range derived from model 3 in Figure 12. The total X- 
ray luminosity is L x = 1.5 x 10 29 ergs s _1 . The X-ray 
is strongest at the head where the shock is the strongest. 
There are strong emission lines from oxygen and iron in 
the spectrum. Current upper limits on the X-ray flux from 
CRL618 are roughly a factor 10 larger than our model flux. 
CRL618 was not detected in a 18697 sec observation with 
the HRI instrument onboard ROSAT (Guerrero, Chu, & 
Gruendl 2000) . If we assume that a minimum of 9 photons 
are required for a 3a detection, we can use the detection 
of BD+30°3639 with HRI (172+/-13 counts in 8740 sec, 
Leahy et al 1998) and its measured X-ray flux 5.7xl0~ 13 
erg cm~ 2 s _1 to put a roughly upper limit on the X-ray 
flux of CRL618 of <1.4xl0~ 14 erg cm" 2 s" 1 (1.7xl0 30 
erg s -1 ). 

The series of optical bow-like structures observed in the 
body of the Wl lobe (as well as other lobes) but not repro- 
duced in the previous models, suggest that the fast wind 
may have time variations in density and velocity. With 
an episodic fast wind, model 4 produces a series of optical 
bow- like structures in the body of the lobe (see Fig. 13). 
The two bright bow- like structures (Bi and B2) at the tip 
trace the newly shocked AGB wind in the outer shell and 
the newly shocked fast wind in the inner shell impacted 
by the internal shocks. However, in contrast to the obser- 
vations that show a series of bow-like structures of similar 
intensity within the lobe, the optical emission of our model 
internal bow-like structures becomes progressively weaker 
away from the source (due to radial expansion and ra- 
diative cooling). The results for the [SII]A6716A/A6730A 
ratio are also inconsistent with the observations of the Wl 
lobe, which show a roughly constant value of 0.5 with fluc- 
tuations from 0.4 to 0.6. In our model, this ratio is (a) 
about 1.5 for the shell along the body of the lobe, (b) 0.6 
in the bow-like structure near the source (and increases 
away from the source). Although there is also emission 
from the inner shell (that is impacted by the wings of the 
internal shocks), this emission does not have the observed 
bow-shaped geometry. 

One way to resolve the above discrepancies between our 
model 4 and the observations is to have a mass-loss rate 
of the fast wind (M/) that decreases with time. How- 
ever, in order to produce adequate optical emission along 
the body of the lobe, M/ at time t = would have to 
be a factor of ~ 10 larger. Another possibility is for the 
underlying driving agent to be an episodic cylindrical jet, 
which has a density constant with the distance. Simula- 
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tions with episodic cylindrical jets have been performed by 
a number of authors in order to reproduce the morphology 
and kinematics of protostellar outflows (Stone & Norman 
1993b; Raga & Cabrit 1993; Biro & Raga 1994; Suttncr 
et al. 1997; Lee et al. 2001). Such jets can also produce 
a series of internal bow-like structures along the body of 
the lobe. Without radial expansion, the emission of these 
structures does not decrease as fast as that in model 4. 
In addition, each of these structures has overall stronger 
emission than its counterpart in our model 4 because it 
closes back to, and is shocked by, the one behind it (see 
Fig. 6 and 7 in Biro & Raga 1994). 

5.2. Kinematics 

The position- velocity (PV) diagrams of optical emission 
for cuts along the axis of the Wl lobe show a nonlinear in- 
crease in velocity with increasing distance from the source 
and a broad range of velocities at the tip with a FWZI up 
to 200 km s _1 (SCSG02). In this section, we present the 
PV diagrams of the [OIJA6300 emission, the most extended 
emission, derived from our simulations and compare them 
to that of the Wl lobe. The diagrams are also degraded 
to low resolutions for comparing with the current observa- 
tions. 

The PV diagrams derived from the simulation of model 
1 are presented in Figure 14. At an inclination i = 30°, the 
newly shocked fast wind is projected at a velocity higher 
than that of the newly shocked AGB wind. At low resolu- 
tion, the diagram shows two elliptical structures attached 
together. These diagrams show a broad range of veloci- 
ties at the tip, qualitatively consistent with that seen in 
the Wl lobe. However, since the emission is mostly from 
the newly shocked AGB wind, the velocity range is about 
v so sin i = 89 km s _1 , much smaller than that seen in the 
Wl lobe. In order to produce the observed velocity range, 
the shock velocity at the tip would have to be 400 km s _1 , 
resulting in a very high temperature for the shell (> 10 6 K) 
As mentioned, the emission from the newly shocked fast 
wind could be as strong as that from the newly shocked 
AGB wind. In that case, the detailed PV structure would 
be similar to that of model 2 (see Fig. 15) and the veloc- 
ity range at the tip would be about Vf s'mi = 150 km s _1 , 
similar to that seen in the Wl lobe. 

The PV diagrams for model 3 are presented in Figure 
16. They show a convex spur structure on the redshifted 
side, i.e., a nonlinear increase in velocity toward the tip 
with increasing distance from the source. There is also 
emission on the blueshifted side. A broad range of ve- 
locities is also seen at the tip. The maximum redshifted 
velocity is v so — ^sinz = 165 km s" 1 while the maxi- 
mum blueshifted velocity is — -^sinz = —55 km s _1 , so 
that the velocity range at the tip is equal to the shock ve- 
locity at the tip independent of the inclination, as found 
in a jet-driven bow shock model by Hartigan, Raymond, 
& Hartmann (1987). As a result, this model produces a 
similar velocity range as that seen in the Wl lobe. 

The PV diagrams for model 4 are presented in Figure 17. 
The PV structures associated with the leading shock are 
similar to that seen in model 1. There are also structures 
associated with the internal shocks projected at a veloc- 
ity higher than that of the leading shock, showing a series 
of increases in velocity. These PV diagrams are different 



from that of the Wl lobe, which shows a single nonlin- 
ear increase in velocity with increasing distance. Hence, 
the bow-like structures seen in the Wl lobe may not be 
associated with internal shocks. 



5.3. Line Ratios and Intensities 

Here we compare the optical emission line fluxes de- 
rived from our simulations to the observations (Table 2). 
Since models 1 and 4 can produce a lobe with a collima- 
tion similar to that of the Wl lobe, we compare these in 
detail to the observations. The fluxes always peak at the 
tip of the lobe in both the simulations and observations. 
However, the ratios of the [OI] to [SII] peak fluxes de- 
rived from our simulations are too low compared to those 
derived from high resolution HST observations, indicat- 
ing that the temperatures of the shocked AGB wind at 
the tip in our simulations (~2x 10 4 K) are always much 
higher than observed (~ 8x 10 3 K) (see Fig. 18). This 
conclusion is further supported by the groundbased obser- 
vations. For example, the ratios of the [OIII] and [Nil] 
total fluxes to the [SII] total fluxes in the tip of the lobe 
in the simulations are significantly higher than observed, 
which also indicates that the temperatures at the tip in 
the simulations are much higher than observed. Notice 
that in the observations, the [OIII] to [SII] ratio indicates 
a significantly higher temperature (<~ 30,000 K) than that 
indicated by the ratios of [Nil] and [OI] to [SII] (~ 15,000 
K). Thus, the [OIII] emission comes from a region hot- 
ter than that of the [Nil] and [OI] emission, consistent 
with an existence of a temperature stratification in the tip. 
Our models also show similar temperature stratification, 
although as mentioned above, the model temperatures are 
higher than observed. 

The high temperatures at the tip of the shell in our 
simulations are produced by the dense and cold jet-like 
structure at the tip as it outruns the shell and shocks the 
AGB wind at the head before the shell does. The shell 
at the tip then interacts with the shocked AGB wind, re- 
sulting in a temperature higher than observed. The jet-like 
structure is formed as the shocked fast wind material flows 
toward the tip (Canto, Tenorio-Tagle, & Rozyczka 1988; 
Frank, Balick, & Livio 1996b). In reality, the cold dense 
jet-like material may mix with the hot shell material at 
the tip of the lobe (e.g., due to instabilities), significantly 
lowering the temperature there. Further simulations with 
high resolution are needed to test this idea. It is interest- 
ing to note that the observed line ratios at the tip of the 
Wl lobe (and thus the inferred temperature) are similar 
to those in the internal shocks in model 4. These internal 
shocks lack the jet-like structure at the tip in our simula- 
tions, thus supporting the idea that the jet-like structure 
is responsible for the high temperatures at the tip. 

Although these models (1 and 4) produce a peak [SII] 
flux comparable to that of the observations, the total [SII] 
fluxes in our simulations are much less than observed be- 
cause of the lack of significant emission within the body 
of the lobe. New simulations will need to address this 
discrepancy. Note that such simulations will also have to 
explain that the observed ratio of the [OI] to [SII] fluxes in 
the body of the lobe is smaller than that at the tip (which 
implies that the body of the lobe is hotter than the tip). 
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5.4. The nature and origin of the collimated fast wind in 

CRL 618 

The fast wind in CRL 618 is very unlikely to be a steady 
wind. This conclusion is most directly supported by the 
presence of the optical bow-like structures within the body 
of the lobe (see discussion in §5.1). Furthermore, CRL 618 
is known to have a fast (<~ 200 km s _1 ), compact (size < 
5") bipolar molecular outflow (Ncri et al. 1992). Only our 
pulsed wind model produces material moving at such high 
speeds close to the central source. 

Since AGB winds are thought to be driven by radiation 
pressure on dust grains (see review by Habing 1996), we 
need to consider the possibility that the fast winds in PPNs 
are also powered by radiation pressure. In our models, the 
momentum flux per steradian of the fast wind is 

^ = %T exp [- 3 ^ 55 P ^M-H^r] (17) 

The luminosity of CRL 618 is L ~ 1O 4 L (Goodrich 
1991), giving a radiation momentum flux per steradian of 
P L = L/4ttc ~ 10 34 erg cm s" 1 yr" 1 sr" 1 . Thus, P f > P L 
within the opening angle of the CFW in our models (see 
Table 1), implying that the CFW can not be driven by ra- 
diation pressure. It is likely that the fast wind is launched 
by magneto-centrifugal forces from a magnetized accretion 
disk and star system (Blackman, Frank, & Welch 2001). 
The accretion disk may be formed in a binary system as 
the material flows from one star to another (Soker & Rap- 
paport 2000). In this case, the gravitational binding en- 
ergy in the accreting material is converted into the kinetic 
energy of the fast wind. 

5.5. Application to other PPN: OH 231.8+4.2 

The simulation of model 2 may provide a clue to the 
formation of another PPN, OH 231.8+4.2. This PPN is 
bipolar, showing two strong and wide bow-like structures 
in Ha at the tip of both the northern and southern lobes 
(Bujarrabal et al. 2002). A fast, collimated, dense CO 
structure is seen inside each lobe along the axis (Sanchez 
Contreras et al. 2000). Recently, it has been proposed that 
the two bow-like structures in both northern and south- 
ern lobes trace the forward and backward shocks driven 
by a dense jet (Bujarrabal et al. 2002). Models with an 
overdense cylindrical jet have difficulty producing as wide 
and strong bow-like emission for the backward shock as 
that for the forward shock due to the fast cooling of the 
backward shock (Raga et al. 1995). In contrast, a fast 
collimated wind can readily produce two strong and wide 
bow-like structures at the tip as those seen in model 2 (see 
Figure 15). Coincidentally, the collimation of the lobe and 
the velocities of both the forward and backward shocks in 
this model are similar to those derived by Bujarrabal et 
al. (2002) for the northern lobe of OH 231.8+4.2. As a 
result, if the two bow-like structures are to be associated 
with the forward and backward shocks, a tenuous colli- 
mated fast wind may be needed. If this is the case, the 
observed dense CO collimated structure would represent 
the dense core of a cool, latitudinally stratified fast wind. 
The lengths of the two lobes are different and may be due 
to an intrinsic difference in the velocity on each side of the 
bipolar fast wind (Soker 2002). 



6. SUMMARY AND CONCLUSIONS 

We have presented a number of simulations of a colli- 
mated fast wind interacting with a spherical AGB wind. 
In our simulations, we have used a dense AGB wind with a 
high mass-loss rate typical of PPNs like CRL 618. The pa- 
rameters of the fast wind are chosen to give a shell velocity 
at the tip of about 150 km s^ 1 , comparable to that seen 
in CRL 618. Under these assumptions, the main results 
from the simulations are the following: 

1. The shell dynamics is determined by the velocity 
of the fast wind. For a 300 km s _1 fast wind, the 
cooling of the shocked fast wind is fast enough to 
remove the heat obtained from the fast wind, so 
that the shell is thin and momentum-driven. For 
a 1000 km s _1 fast wind, however, the shocked 
fast wind does not cool fast enough and expands 
sideways into a hot cocoon surrounding the fast 
wind. The fast wind itself becomes a ripple jet, 
so that the shell is essentially driven by a jet and 
its kinematics can be described by a ballistic bow 
shock model. 

2. Although the AGB wind is spherical, the shell 
driven by a collimated fast wind has the shape of 
a highly collimated lobe. For a shell driven by 

a 300 km s _1 fast wind, the collimation of the 
shell increases roughly linearly with the decreasing 
opening angle of the fast wind. However, for a shell 
driven by a 1000 km s _1 fast wind, the collimation 
increases slowly with the decreasing opening angle. 

3. A time- varying velocity of the fast wind produces 
a series of internal shock pairs interacting with the 
inner surface of the shell. Due to radial expansion 
of the fast wind, the density of the internal shocks 
decreases rapidly with distance. 

We have also derived various emission diagnostics from 
our simulations. The main results are the following: 

1. For a 300 km s _1 fast wind, the shell is composed 
of the shocked AGB wind (outer shell) and the 
shocked fast wind (inner shell). The optical 
emission arises mostly from the newly shocked 
AGB wind at the tip of the lobe, forming a 
bow-like structure. If the velocities of the shocked 
fast wind and the shocked AGB wind are similar, 
two bow-like structures are seen at the tip because 
the newly shocked fast wind also contributes 
significantly to the emission. 

2. For a 1000 km s _1 fast wind, the shell is 
dominantly composed of the shocked AGB wind, 
and the optical emission again forms a bow-like 
structure at the tip. The shocked fast wind, 
however, forms a tenuous and hot, X-ray emitting 
structure interior to the shell. 

3. When the velocity of the fast wind is time-varying, 
bow-shaped optical emission structures are seen 
associated with the internal shocks along the body 
of the lobe. This emission decreases rapidly away 
from the source. 
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4. The position-velocity (PV) diagrams from all our 
models show a broad range of velocities at the tip. 
However, the velocity range at the tip depends on 
the shell dynamics and the relative contributions of 
the shocked fast wind and the shocked AGB wind. 
For a momentum-driven shell, the velocity range 
depends on the inclination. For a shell driven by a 
ballistic bow shock, however, the velocity range is 
independent of the inclination and is equal to the 
shock velocity at the tip. When the velocity of the 
fast wind is time-variable, a series of increases in 
the velocity are seen in the PV diagrams. 

Comparing our simulations to the observations of the 
PPN CRL 618, we find that 

1. A 300 km s _1 collimated fast wind with an 
opening angle of 10° can readily produce a highly 
collimated lobe similar to the Wl lobe in this 
object. However, a 1000 km s _1 fast wind, even 
with an opening of 5°, has difficulty producing 
such a highly collimated lobe. 

2. Our models can produce an optical bow-like 
structure at the tip similar to that seen in the 
Wl lobe. However, they have difficulty producing 
several additional bright optical bow-like structures 
seen in the body of the lobe. A fast wind with a 
less steep (than r~ 2 ) radial density gradient may 
help to remove this difficulty. 

3. The PV diagrams derived from our simulations 
are qualitatively consistent with the observations. 



However, the observed velocity range at the tip of 
the Wl lobe is larger than the momentum-driven 
shell model and may favor a shell driven by a 
ballistic bow shock. 

4. The [SII]A6716A/A6730A ratios at the tip of the 
lobe in our simulations are all about 0.5 (implying 
an electron density of ~ 10 cm -3 ), as observed. 

5. Specific line ratios indicate that the temperatures 
at the tip of the lobe in our simulations are higher 
than observed. This may result from insufficient 
mixing of hot and cold material in the tip in 
our simulations. A comparison of the inferred 
temperatures from the [OIII]/[SII] ratio to that 
from the [NII]/[SII] and [OI]/[SII] ratios indicates 
a temperature stratification in the tip, both for the 
simulations and observations. 

6. The collimated fast wind is unlikely to be steady 
and is not radiatively driven. 
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APPENDIX 

EMISSION CALCULATION 

The optical emission is assumed to be a thermal equilibrium (LTE) optically thin emission, so that the emissivity is 
given by 

xfif u nA u iE u i 

6 = ^ (A1) 

where x, /„, A u i, E u i are the abundance of the element relative to hydrogen nuclei, ionization fraction of the ions, 
fraction of the ions in the upper excited state, Einstein A coefficient, and energy separation between the upper and lower 
excited states, respectively. The abundances relative to hydrogen nuclei are assumed to be 4.17xl0~ 4 , 9.77xl0~ 6 and 
1.28xl0 -4 for O, S and N, respectively. The ionization fractions of [OI] and [OIII] as a function of temperature are taken 
from Nahar (1999), [Nil] from Nahar & Pradhan (1997) and [SII] from Shull & van Steenberg (1982). /„ is obtained by 
solving the level populations with the relevant atomic quantities taken from Mendoza (1983). The X-ray emission is also 
assumed to be an equilibrium optically thin emission. It is calculated with the Raymond & Smith (1977) model and a 
solar abundance from Anders & Grevesse (1989). 



11 



REFERENCES 



Anders, E. & Grevesse, N. 1989, Geochim. Cosmochim. Acta, 53, 
197 

Balick, B. & Frank, A. 2002, ARA&A, 40, 439 
Balick, B., Preston, H. L., & Ickc, V. 1987, AJ, 94, 1641. 
Balick, B., Rugers, M., Terzian, Y., & Chengalur, J. N. 1993, ApJ, 
411, 778. 

Biro, S. and Raga, A. C. 1994, ApJ, 434, 221 
Blackman, E. G., Frank, A., & Welch, C. 2001, ApJ, 546, 288 
Blandford, R. D. & Payne, D. G. 1982, MNRAS, 199, 883 
Blondin, J. M. & Lundqvist, P. 1993, ApJ, 405, 337. 
Blondin, J. M., Fryxell, B. A. & Konigl, A. 1990, ApJ, 360, 370 
Borkowski, K. J., Blondin, J. M., & Harrington, J. P. 1997, ApJ, 
482, L97. 

Bujarrabal, V. & Alcolea, J. 1991, A&A, 251, 536 
Bujarrabal, V., Castro-Carrizo, A., Alcolea, J., & Sanchez 

Contreras, C. 2001, A&A, 377, 868 
Bujarrabal, V., Alcolea, J., Sanchez Contreras, C, & Sahai, R. 

2002, A&A, 389, 271 
Canto, J., Tenorio-Tagle, G., & Rozyczka, M. 1988, A&A, 192, 287. 
Chevalier, R. A. & Luo, D. 1994, ApJ, 421, 225 
Chu, Y., Jacoby, G. H., & Arendt, R. 1987, ApJS, 64, 529 
Cliffe, J. A., Frank, A., Livio, M., & Jones, T. W. 1995, ApJ, 447, 

L49 

Corradi, R. L. M. & Schwarz, H. E. 1995, A&A, 293, 871 
Dalgarno, A. & McCray, R. A. 1972, ARA&A, 10, 375 MNRAS, 
314, 241 

Dwarkadas, V. V. & Balick, B. 1998, ApJ, 497, 267. 
Dwarkadas, V. V., Chevalier, R. A., & Blondin, J. M. 1996, ApJ, 
457, 773. 

Dyson, J. E. & Williams, D. A. 1980, The Physics of the Interstellar 

Medium (New York: Wiley) 
Frank, A. 1999, New Astronomy Review, 43, 31 
Frank, A. & Mcllcma, G. 1994, ApJ, 430, 800. 
Frank, A. & Mellema, G. 1996, ApJ, 472, 684 
Frank, A., Balick, B., & Livio, M. 1996, ApJ, 471, L53. 
Frank, A., Ryu, D., & Davidson, K. 1998, ApJ, 500, 291 
Garcia-Segura, G. 1997, ApJ, 489, L189 
Gardiner, T. A. & Frank, A. 2001, ApJ, 557, 250 
Goodrich, R. W. 1991, ApJ, 376, 654 

Guerrero, M. ;., Chu, Y., & Gruendl, R. A. 2000, ApJS, 129, 295 
Habing, H. J. 1996, A&A Rev., 7, 97 

Hartigan, P., Raymond, J., & Hartmann, L. 1987, ApJ, 316, 323 

Hartmann, L. & Kenyon, S. J. 1996, ARA&A, 34, 207 889 

Icke, V., Balick, B., & Frank, A. 1992, A&A, 253, 224. 

Kahn, F. D. 1983, IAU Symp. 103: Planetary Nebulae, 103, 305. 

Knapp, G. R. & Morris, M. 1985, ApJ, 292, 640 

Konigl, A. & Pudritz, R. E. 2000, Protostars and Planets IV, 759 



Kwok, S., Purton, C. R., & Fitzgerald, P. M. 1978, ApJ, 219, L125. 
Kwok, S., Su, K. Y. L., & Hrivnak, B. J. 1998, ApJ, 501, L117. 
Lee, C, Stone, J. M., Ostriker, E. C, & Mundy, L. G. 2001, ApJ, 
557, 429. 

MacDonald, J. & Bailey, M. E. 1981, MNRAS, 197, 995 
Mellema, G. 1995, MNRAS, 277, 173. 

Mcndoza, C. 1983, IAU Symp. 103: Planetary Nebulae, 103, 143. 

Nahar, S. N. 1999, ApJS, 120, 131. 

Nahar, S. N. & Pradhan, A. K. 1997, ApJS, 111, 339. 

Neri, R., Garcia-Burillo, S., Guelin, M., Cernicharo, J., Guilloteau, 

S., & Lucas, R. 1992, A&A, 262, 544 
Neri, R., Kahanc, C, Lucas, R., Bujarrabal, V., & Loup, C. 1998, 

A&AS, 130, 1. 

Ostriker, E. C, Lee, C, Stone, J. M., & Mundy, L. G. 2001, ApJ, 
557, 443. 

Raga, A. C. & Cabrit, S. 1993, A&A, 278, 267 
Raga, A. C, Taylor, S. D., Cabrit, S., & Biro, S. 1995, A&A, 296, 
833 

Raymond, J. C. & Smith, B. W. 1977, ApJS, 35, 419 
Rozyczka, M. & Franco, J. 1996, ApJ, 469, L127 
Sahai, R. 2001, Post-AGB Objects as a Phase of Stellar Evolution, 
53. 

Sahai, R. 2003, IAU, in press 

Sahai, R. & Trauger, J. T. 1998, AJ, 116, 1357 

Sahai, R. et al. 1998, ApJ, 493, 301. 

Sanchez Contreras, C. & Sahai, R. 2001, ApJ, 553, L173. 
Sanchez Contreras, C, & Sahai, R. & Gil De Paz, A 2002, ApJ, in 
press 

Sanchez Contreras, C, Bujarrabal, V., Miranda, L. F., & 
Fernandcz-Figueroa, M. J. 2000, A&A, 355, 1103 

Schwarz, H. E., Aspin, C, Corradi, R. L. M., & Reipurth, B. 1997, 
A&A, 319, 267 

Shu, F., Najita, J., Ostriker, E., Wilkin, F., Ruden, S., & Lizano, S. 

1994, ApJ, 429, 781 
Shull, J. M. & van Steenberg, M. 1982, ApJS, 48, 95. 
Soker, N. 2002, ApJ, 568, 726 
Sokcr, N. & Rappaport, S. 2000, ApJ, 538, 241 
Steffen, W. & Lopez, J. A. 1998, ApJ, 508, 696 
Stone, J. M. & Norman, M. L. 1993, ApJ, 413, 210 
Stone, J. M. & Norman, M. L. 1994, ApJ, 420, 237 
Su, K. Y. L., Volk, K., Kwok, S., & Hrivnak, B. J. 1998, ApJ, 508, 

744. 

Suttner, G., Smith, M. D., Yorke, H. W. & Zinncckcr, H. 1997, 
A&A, 318, 595 

Trammell, S. R. 2000, ASP Conf. Ser. 199: Asymmetrical Planetary 
Nebulae II: From Origins to Microstructures, 147. 



12 



Table 1 
Parameters For Fast Winds 



Model 


Mf 

(10" 7 M yr- 1 ) 


v fo 

(100 km s- 1 ) 


°f 
(degree) 


R b 

(io- 3 ) 


K 

(%) 


Pfo/PL 


1 


25 


3 


10 


7.6 


86.5 


158 


2 


100 


3 


20 


30.2 


39.6 


158 


3 


2.5 


10 


10 


7.6 


86.5 


53 


4* 


25 


3 


10 


7.6 


86.5 


158 



a Opening angle 

b Total mass- loss rate within 10° from the pole/Total mass- loss rate 

c Momentum flux per steradian of the fast wind at the pole/Radiation momentum 
flux per steradian 

* Model 4 is Model 1 with time variation in the density and velocity. 



Table 2 
Emission line diagnostics 



Peak Intensity 




Model 




Observations" 




1 2 3 


4 




[OI]A6300A C 


0.5 0.9 1.3 


0.53 


5.4 


[SII]AA6730,6716A b 


l.le-13 6.1c-13 l.le-12 


8.1c-14 


1.4c-13 


[NII]A6583A C 


5.0 6.8 9.3 


3.6 




[OIII]A5007A C 


5.6 7.5 19.7 


90 






Total Flux in the Peak d 








Model 




Observations 6 




1 2 3 


4 




[OI]A6300Af 


0.6 0.8 0.8 


0.7 


1.5 


[SII]AA6730,6716A^ 


6.6e-15 3.2e-13 3.7c-13 


5.1c-15 


l.le-14 


[NII]A6583A 9 


3.9 5.0 6.5 


2.0 


1.3 


[OIII]A5007Af 


2.1 3.2 9.8 


0.7 


0.07 


a Adopted from the HST images in Sahai et al. (2003) 




b Pcak Intensity in 


unit of erg s 1 cm 2 arcsec 


-2 





Teak Intensity normalized to that of [SII]AA6730,6716A 
d Mcasured in a 2" x2" box covering the tip of the Wl lobe. 



e From groundbased observations (SCSG02 and priv. comm.) 

f Total flux in the peak in unit of erg s~ x cm~ 2 . 

STotal flux in the peak normalized to that of [SII]AA6730,6716A 
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Fig. 1. — Simulation of model 1 as the fast wind propagates to a distance of ~ 10 17 cm, a length similar to that of the lobes in CRL 
618. The age of the fast wind is 177 years. In this model, Mj = 2.5 X 10 — 6 Mq yr~ 1 , vj = 300 km s~ 1 , and dj = 10°. (a) shows the 
contours of the fast wind tracer, c, superposed on the density distribution. Contours are 0.1, 0.5 and 0.9. The dashed line is the shell structure 
derived from the momentum-driven shell model with Equation 10. (6) and (c) show the contours of c w and c a , respectively, superposed on 
the density distribution. c a and c w trace respectively the flows of the AGB wind and the fast wind emanating within the angular segment 
5°< 8 <7.5°. Contours are 0.1, 0.5 and 0.9. (d) shows the velocity distribution superposed on the pressure distribution, (e) shows the 
temperature distribution. The gray-scale wedges along the right side indicate the values for the gray-scale images. (A), (B) and (C) show 
respectively the blow-ups for the regions in the boxes A, B and C. Regions i, ii, iii, iv and v indicate unshocked fast wind, shocked fast wind, 
mixture of shocked fast and AGB winds, shocked AGB wind and unshocked AGB wind, respectively. Contact discontinuity is located in 
region iii. 



Fig. 2. — Transverse velocity (vr) and longitudinal velocity (v z ) of the shell material for model 1 at 177 years. The gray-scale indicates 
the mass column density on a linear stretch. The solid lines indicate the R and z components of the mean shell velocity obtained from the 
simulation. The dashed lines indicate the R and z components of the velocity of the newly shocked fast wind derived from Equation 11. 



Fig. 3. — Coordinate system used for the shell dynamics. The thick line indicates the shell. Vf, v s and r a are the velocity of the fast wind, 
the shell velocity and the shell radius, respectively, ct is the angle between the shell normal and the shell velocity. 



Fig. 4.— Simulation of model 2. The age is 188 years. In this model, M f = 1 x 10~ o M© yr -1 , v fo = 300 km s" 1 , and 6 f = 20°. The 
captions are the same as Figure 1. 



Fig. 5. — Simulation of model 3. The age is 180 years. In this model, Mf = 2.5 X 10~ 7 M yr _1 , v fo = 1000 km s _1 , and 6/ = 10°. The 
captions are the same as Figure 1, except that the dashed line in (a) is the shell structure derived from the ballistic bow shock model given 
by Ostriker et al. (2001). 



Fig. 6. — Transverse velocity (vr) and longitudinal velocity (v z ) of the shell material for model 3 at 180 years. As in Fig. 2, the gray-scale 
indicates the mass column density on a linear stretch. The solid lines indicate the R and z components of the mean shell velocity, while the 
dashed lines indicate the R and z components of the immediate postshock velocity of the AGB wind derived from the ballistic bow shock 
model. 



Fig. 7. — Simulation of model 4. The age is 166 years. In this model, the fast wind has the same parameters as that in model 1 but has a 
time-variation with A = 0.5 and P = 22 years (see Equation 5). The captions are the same as Figure 1, except that the dashed lines in (a) 
are derived from Equation 13. 



Fig. 8. — HST images of CRL 618 in [OI] and [SII] emission (logarithmic stretch). The cross indicates the approximate position of the 
central source. The gray-scale wedges along the right side indicate the values for the gray-scale images. The intensity of the image has a unit 
of erg s _1 cm -2 arcsec 2 . 



Fig. 9. — Integrated optical emission and line ratio of [SIIJA6716A/A6730A for model 1 at an inclination of i=30°. The optical emission is 
in logarithmic scale with an angular resolution similar to that of HST images. The gray-scale wedges along the right side indicate the values 
for the gray-scale images. 



Fig. 10. — Same as that in Figure 9 except for model 2. 



Fig. 11 — 



Same as that in Figure 9 except for model 3. 
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Fig. 12. — Integrated X-ray emission and spectrum for model 3 in the 0.2-1.5 keV energy range at an angular resolution of 0.5" and an 
energy resolution of 50 eV. The gray-scale wedge along the right side indicates the values for the gray-scale image. In the spectrum, there are 
strong emission lines from oxygen at about 570 and 653 eV, and from iron at about 730 and 826 eV. 



Fig. 13. — Same as that in Figure 9 except for model 4. 



Fig. 14. — Position- velocity diagrams of the [OI] emission for model 1 at an inclination of i=30°. PV diagram on the left has an angular 
resolution of 0.05" and a velocity resolution of 10 km s _1 . PV diagram on the right is degraded to an angular resolution of 1" and a velocity 
resolution of 50 km s _1 for comparing with current groundbased observations. 



Fig. 15. — Same as that in Figure 14 except for model 2. 



Fig. 16. — Same as that in Figure 14 except for model 3. 



Fig. 17. — Same as that in Figure 14 except for model 4. 



Fig. 18.— Line ratios of the [OI]A6300A, [NII]A6583A and [OIII]A5007A emission to the [SII]AA6730,6716A emission. Solid lines are for 
number density=2xl0 3 cm~ 3 , dashed lines for 10 4 cm -3 , and dotted-dashed lines for 5xl0 4 cm -3 . The line ratios are derived assuming 
equilibrium optically thin emission (see Appendix). 
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